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The problem of two electrons in a square billiard interacting via a finite-range repulsive Yukawa 
potential and subjected to a constant magnetic field is considered. We compute the energy spectrum 
for both singlet and triplet states, and for all symmetry classes, as a function of the strength and 
range of the interaction and of the magnetic field. We show that the short-range nature of the 
potential suppresses the formation of "Wigner molecule" states for the ground state, even in the 
strong interaction limit. The magnetic susceptibility x(B) shows low-temperature paramagnetic 
peaks due to exchange induced singlet-triplet oscillations. The position, number and intensity of 
these peaks depend on the range and strength of the interaction. The contribution of the interaction 
to the susceptibility displays paramagnetic and diamagnetic phases as a function of T. 
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I. INTRODUCTION 

The study of mesoscopic systems has proved to be a 
rich field to investigate explicit manifcstations_of quan- 
tum properties in nanometer and micron scalestl. In such 
systems, the electron coherence length scales and mean- 
free paths are in general larger than the typical sample 
sizes, so that the underlying classical electronic motion 
plays an important role. The nature of the classical mo- 
tion, regular, mixed or chaotic, reflects itself on some 
of the quantum properties of the system, particularly in 
the energy level distribution. These features have long 
been studied for non-interactinga and weakly interacting 
systems!! 

In quantum dots, where few electrons are laterally 
confinedta, the electron-electron interaction is usually 
very efficiently shielded by the positively charged fixed 
ions and other effects, so that the independent electron 
gas theory can often be. .-use d j to , understand the basic 
features of the systemHBfl&BliSLil. However, the resid- 
ual interaction that survives the shielding can sometimes 
play important roles. In recent years, much attention has 
been giy nn .. tQ u j itpjapf jon- induced effects in mesoscopic 
systemsBE3tJ£aO'E2lEil. These effects are particularly 
important in the large-dot regime, where the electron- 
electron Coulomb interaction overcomes the kinetic en- 
ergy, forcing the ground-state into a Wigner molecule 
type of configuratior£3. 

In mesoscopic systems the electronic interaction is usu- 
ally not well approximated by a bare long-range Coulomb 
force, due exactly to screening effectsEj. The strength 
and range of the residual interaction, or the efficiency 
of the screening, depend on many parameters, like the 
electron density and size of the dot. It is therefore im- 
portant to understand the effects of the interaction as a 
function of its effective intensity and range. In this work 
we give a contribution in this direction presenting exact 



results for the problem of two electrons in a square quan- 
tum dot interacting via a repulsive finite-range Yukawa 
type of interaction, V(r) = V~oe~ ar /r, and subjected to 
a uniform and constant magnetic field of strength B ap- 
plied perpendicular to the dot. This model system was 
inspired by the experimental work of Levy et al3 where 
the orbital magnetic susceptibility was measured for an 
ensemble of square dots containing of the order of a thou- 
sand electrons each. Although the independent particle 

gpmirlassiral thpnry pyplainpH most of the experimen- 
tal findingsElM'BBy, the behavior of the susceptibil- 
ity with the temperature does not come out correctly in 
this approach and still puzzles the theorists. The idea 
that the slow decay of the susceptibility observed exper- 
imentally (as opposed to the exponential decay expected 
from the semiclassical theory), could be duetto electron- 
electron interaction was first investigated inEI for a weak 
contact (Dirac delta) type of interaction using pertur- 
bation methods. In this article we study the effects of 
electronic interaction a much simpler system, with only 
two electrons, but we present exact (numerical) results. 

The choice of an Yukawa type of potential allows us 
to interpolate between the pure Coulomb (a ~ 0) and 
short range interactions. Besides, the calculation of the 
Hamiltonian matrix elements can be reduced to one- 
dimensional integrals, which can be calculated numeri- 
cally This allows us to compute the energy spectrum 
for the four rotational symmetry classes as a function of 
the interaction strength Vq and range 1/a for both sin- 
glet and triplet states. We also consider these results 
as a function of a constant magnetic field of strength B 
applied perpendicular to the square. We compute the 
magnetic susceptibility at finite values of the magnetic 
field and temperature via the partition function. 

Our main results can be summarized as follows: (1) 
Vq introduces avoided crossings between the energy lev- 
els within each symmetry class, one of the signatures 
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of quantum chaos; (2) a has a very important role in 
determining the probability profile of the ground state, 
suppressing in some cases the Wigner molecule type of 
behavior even for strong interactions; (3) The effect of 
the interaction on the magnetic susceptibility x(B) de- 
pends on a. In particular, for large magnetic fields, 
singlet-triplet oscillations of the ground state level lead 
to paramagnetic fluctuations on the two-electron suscep- 
tibility x(-B), in contrast to the non-interacting diamag- 
netic susceptibility. The position and intensity of these 
peaks change with the range of the interaction. ; (4) The 
contribution to the susceptibility induced by interaction 
at zero magnetic field shows paramagnetic and diamag- 
netic phases as a function of the temperature. This type 
of behaviorr-,has also been found for weak Dirac delta 
interactions!! 

This paper is organized as follows: in section II we 
describe the system in detail. We discuss its symme- 
try properties and compute the matrix elements of the 
Hamiltonian. In section III we present numerical results 
for the energy spectrum and for ground-state electronic 
density as a function of the strength and range of the in- 
teraction. In section IV we consider the magnetic proper- 
ties of the system and in section V we discuss our results. 



II. HAMILTONIAN AND MATRIX ELEMENTS 

We consider a system where two electrons are confined 
in a square-shaped two dimensional billiard of size L in- 
teracting via an Yukawa type of potential and subjected 
to a uniform and constant magnetic field of strength B 
applied perpendicular to the dot. The Hamiltonian is 
given by 

H = TTlE (P*< + eB y*l 2 ) 2 + (Ptf ~ eBxt/2) 2 ) 



i=l,2 



e -a|ri-r 2 | 

+ V wa lls + V — — (1) 

\n - r 2 \ 

where m* is the quasiparticle electron mass and V wa ii s is 
the square well potential. The eigenfunctions of a single 
particle in this square dot with zero magnetic field are 
given by a normalized product of sine functions: 



fmn{x,y) 
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and the eigenenergies are simply 
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The square billiard is a highly symmetric system. It is 
invariant under the action of 8 symmetry operations (4 
rotations plus 4 reflections) which form the C^ v symmetry 
group. When the time-reversal symmetry is broken (e.g. 
by the application of a magnetic field), the Hamiltonian 
is no longer invariant under reflections. The group then 



reduces to C4, formed by the four rotations generated 
by C4 (rotation by -| ). The symmetric eigenfunctions 
can be written as a linear combination of a particular 
eigenfunction and its symmetry-related counterparts: 



i>(x, V) = ^{x, y) + Ci(p(x, y) + C%ip(x, y) + C%(p(x, y) . 

(4) 

Rotating tp leads to = e* e ip with (e lf >) 4 = 1. This, 
in turn, leads to four solutions for e l6 , namely +1, — 1, +i 
and —i. We can thus separate the general eigenfunc- 
tions ((2J) in four "classes," (or representations) using the 
group's character tableE3, as follows: 



(x, y) if n = m (both odd) 

^= (ip mn (x,y)(±)tp nm (x,y)) 

if n 7^ m and both odd (even) 

(x,y) if n = m (both even) 

((Pmn{x,y)(T)ip„ m (x,y)) 

if n 7^ m and both odd (even) 

,X+i)( \ _ f 7f (<Pmn(x,y) ± i<Pnm(x,y)) 
Wmn (X,y) — < 

it m even (odd) n odd (even) 
.(-<)/_ , A _ J 71 {Vmn{x,y) T i<Pnm(x,y)) 

(_ 11 m even (odd) n odd (even) 

(5) 

This equation can be written in a more compact form 

as 



(G) 



where I = (m,n) and I = (n, m). is either 1 or 

(c) 

and Si is 0, ±1 or ±i, depending on the symmetry class 
(C) and on I (whether (m,n) is odd/even and whether 
m = n or m 7^ n). 

Finally the two-particle orbital eigenfunctions are sym- 
metrized (S) or anti-symmetrized (A) combinations of 
one-particle orbital eigenfunctions: 



(S,A),(C 1 ®C 2 ) 



\ (^ Cl) (rl)< C2) (r- 2 ) ± ^\f^\n. 



V2 



(7) 



The orbital eigenfunction is symmetrized if the electrons 
are in the singlet spin state and anti-symmetrized if they 
are in the triplet spin state. 

The symmetry group of the two-particle system is C4 £§> 
C4 and the eigenfunctions still separate in four symmetry 
classes. The two-particle (2p) class is defined by the total 
phase gained under the action of an element of the C4 (£> 
C4 group (E C4, C48 C4 and so on). This phase is 
simply the product of the one-particle (lp) phases in the 
representation shown in (||) . The 2p class is thus obtained 
in a simple manner by "multiplying" the lp classes. For 
instance, two lp states of class (—1) form a 2p state of 
class (+1) [pictorically, (+1) "=" The same 
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happens with a Ip state of class combined with other 
from class (— i). On the other hand, two lp (+£) states 
form a (—1) 2p state and so on. 



action term dominates over the kinetic term for large L. 
Thus, we define an "effective potential strength" V$ that 
grows linearly with the dot size L: 



A. The screened Coulomb interaction 

For the electron-electron interaction we have used an 
''Yukawa type" short-range screened Coulomb potential 



a\fi— ri\ 



47re e r In - fil 



(8) 



where 1/a is the interaction range and e r is the dielectric 
constant of the two-dimensional electron gas (in case of 
a GaAs 2DEG, e r = 10.9). The reason for this particular 
choice of screening is two-fold. First it interpolates be- 
tween the pure Coulomb case and localized interactions. 
Also it gives an effective "interaction length" (a -1 ) which 
is easy to control. Second, the 1/r dependence is main- 
tained with the screening appearing as an exponential (as 
opposed to a power of This facilitates enormously 

the calculation of the matrix elements, as we show in the 
appendix. The range a will be considered here as a free 
parameter. 

The Hamiltonian for the two electrons without the 
magnetic field is given by 



H 



i=l 



+ v(fi,f2) . 



(9) 



Since the kinetic energy (0) scales with 1/L 2 and the 
interaction scales with 1/L, the electron-electron inter- 



ne) = L 
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so that we can write the matrix elements of H in the 
non-interacting basis in units of (ir 2 h 2 ) / (2m* L 2 ): 



<42\m2 



TY 2 ti 

2m 
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(11) 



where = (mf + nf + m 2 , + n 2 ,) is the kinetic en- 

ergy in units of (7r 2 h 2 )/(2m* L 2 ). This defines another 
free parameter, Vq (or, equivalently, L) that controls the 
relative interaction strength. 

The next step is to calculate the matrix elements of the 
potential in the two-particle eigenfunction basis defined 
by Eq.(^). The repulsive potential does not break the 
Ca <8> Ca rotational symmetry of the Hamiltonian, since 
it depends only on the distance between the electrons. 
Therefore, the interaction matrix is block-diagonal in this 
representation, i.e. the matrix elements 



V, 



_ / .S, A, (Class) 



L l l 2 



S,A, (Class) 



(12) 



are non-zero only inside the same symmetry block. 
For totally symmetric (antisymmetric) eigenfunctions 
^hhi' i' b rea ks into a sum (difference) of the a direct and 

an exchange term. The general expression for the V ma- 
triz elements is given in terms of the general one-particle 
states Eq.(0) as 



F w F w F p F w x 



\i + s[ A) s ( 2 B) s[ c) s { 2 D) ) < hh\V\l[l' 2 > +(s[ A) *s { 2 B) *s[ c) * + s ( 2 D) ) < hh\V\l'A > 
Us{ A) *s { 2 B) *s {D) * + sp) < hhiviU > +(s[ A) *sp*s ( f >* + 4 S) ) < hh\V\l'il2 > 

KS (B)* s (c)* s tn )* + 5 U) } < - lM y^ 2 > +i s[ A) *sl B > + S^SW) < hh\V\?A > 



-(s[ A) * s { ? } * 



+ s^s^)<i 1 i 2 \v\i' i r 2 > + (s ( 2 B) *s\?> 



V 2 

+ s (A) s { f ) )<hi 2 \v\i' 1 T 2 > 



The terms < hhlVl^^ > can be written explicitly with help of Eq.(||) as 



(13) 



16 f . -Km\X\ . ivmyi . 7rm 2 x 2 . 7rn 2 j/2 , rn ^ . imi^xr . im^i . irm 2 x 2 . 7rn 2 t/2 ,2 ^ ,2 ^ 
yj / sm sm — - — sm sm — - — V(\r2— ri\) sm sm — - — sm sm — - — a rid r 2 

1J I Li Li Li Li Li Li Li Li 



(14) 



r 



Eq.(|13[) can be further simplified using the property 



(hh\V 



(hh\V 



1^2 j. The integrals in Eq.(|lJ) 

can be evaluated by switching to relative polar (r, 9) and 
center-of-mass coordinates. Thanks to the exponential 



form of the screening, three of the four integrals in (|1J) 
can be done analytically. The remaining integral, over 
the relative polar angle 9, is performed numerically. Most 
of the direct and exchange elements involve less than 16 
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integrals. The number is actually where N is the 
number of one-particle states with m = n involved in 
either one of the two-particle functions. All these facts 
reduce the number of numerical integrals to be evaluated. 
The details of this calculation are on Appendix A. 

B. The magnetic field 



These are typical of systems with GOE-type level spac- 
ing distribution. Although the number of levels does not 
allow for a precise statistical analysis of the spectrum, the 
level-spacing histograms (not shown) display a distinctive 
difference between the non-interacting (Poisson-like) and 
the interacting (GOE-like) cases. As the relative range 
decreases (aL increases) the levels become less sensitive 
to Vq and the avoided crossings narrow. 



For B 0, there are additional terms in the kinetic 
matrix element of ([!]) proportional to B and B 2 . These 
terms lead to integrals combining sine and cosine func- 
tions and powers of x and y, which can all be done ana- 
lytically. 

The terms linear in B (involving p y and p y ) contribute 
imaginary parts to the matrix elements, breaking the de- 
generacy of the and (— i) symmetry classes. This is 
a consequence of the time-reversal symmetry breaking. 



III. EFFECTS OF INTERACTION: STRENGTH 
AND RANGE 

In this section, we show the numerical results obtained 
with the exact diagonalization of the two-particle inter- 
acting Hamiltonian without the magnetic field. We dis- 
cuss the effects of the two independent parameters of our 
model: the strength Vq and the range a. To change the 
intensity of the interaction relative to the kinetic energy 
we need to change L (see Eq.(10)). However, changing L 
changes the energy levels even if Vq = 0. Therefore, in or- 
der to focus on the changes induced only by the potential, 
we shall measure the energy in units of 2 m'L' 2 throughout 
this section. In these units, the non-interacting eigen- 
energies are independent of L\ the ground state energy, 
in particular, is equal to 4. 



A. Energy levels 

We first consider the effects of the interaction strength 
Vq. To increase Vq relative to the kinetic energy we need 
to increase L. That, however, decreases the effective 
range of the interaction. To keep the ratio between range 
and size of the dot fixed, and concentrate on the effects 
of the potential strength, we shall keep the relative range 
l/(aL) fixed as L (or Vb) is changed. When l/(aL) > 1, 
the electrons "feel" the presence of each other all over 
the dot. For l/(aL) < 1, the interaction is more local- 
ized and the interaction range reduced. 

Figure [l] shows the two-particle energy levels as a func- 
tion of the interaction strength for different values of aL. 
All energy levels are shown, for the four symmetry classes 
of both singlet and triplet configurations. The first panel 
shows the Coulomb case-- a — 0. The interaction in- 
duces singlet-triplet gapstZl and removes several degen- 
eracies in the energy levels. It also promotes level repul- 
sion ( "avoided crossings" ) within each symmetry class. 



B. Ground-state properties 

Recent works have investigated the formation of 
"Wigner molecule" type of ground states .-iii-.palygo- 
nal quantum-dots in the low density limitllSEJO. In 
this limit, the Coulomb interaction between the elec- 
trons dominates over the kinetic energy (the so called 
"large r s " limit) and the ground state electron density 
shows pronounced peaks near the corners of the polygo- 
nal boundaryo. 

We have addressed the question of whether the finite- 
range character of the repulsive potential would change 
such configuration. The low-density limit can be ap- 
proached by making Vq — > oo. However, as dis- 
cussed above, as the dot size L increases, the interaction 
strength increases but the effective interaction range a -1 
decreases. Figure ^| shows that, depending on the value 
of a -1 , the Wigner molecule state can be suppressed, 
even for large dots. This figure shows the ground-state 
electronic density 



J |*o( 



(r!,r~5)| df 2 



(15) 



for L = 10, 100 and 1000 nm and a" 1 = 10, 100 and 
1000 nm. Each column represents dots with the same 
width L but with different values of a. Each line has the 
same value of a but different sizes. The product aL is 
constant along the diagonal. Even for the largest dot, 
with lOOOnm, the Wigner molecule state is completely 
suppressed for aL < 10 (first two plots on the last col- 
umn). Only when aL = 1 (last plot on the last column) 
does the electron density show pronunced peaks near the 
corners. 

Figure || shows two examples of states with peaks near 
the corners in the strong interaction limit. We see that, 
although the electronic density of the two states looks 
similar, the spatial correlations are very different, reflect- 
ing the fact that one of them is a singlet and the other 
a triplet. Fixing ri at the center of one of the peaks, 
Fi = (0.2,0.2), the probability density |*(ri = Fi,r 2 )| 2 
shows two peaks along the diagonal for the singlet state 
and only one peak on the opposite corner for the triplet 
state, since -0 = for ri = r 2 in this case. The two- 
particle configuration is shown schematically for both 
cases. 
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IV. EFFECTS OF THE INTERACTION IN THE 
MAGNETIC PROPERTIES 

A. Energy levels 

Figure [| shows the first energy levels as a function of 
the magnetic field for L — 200 nm and different values of 
a. The first plot shows the non-interacting case (Vq = 0), 
where the singlet and triplet two-particle levels are de- 
generate. Notice that the symmetry classes and (— i) 
are no longer degenerate for B / 0. Also, when the elec- 
tronic interaction is switched on (Figs.^(b) and (c)), the 
singlet-triplet degeneracy is broken. The combination of 
these two effects leads to singlet-triplet crossings in the 
ground state for magnetic fields of the order of a few 
Tesla. This kind of oscillation has been studi 



ously both theoretically (for the Coulomb case)li 
and experimentallyESE3. 

The role of the potential range a can also be seen from 
these figures. In Fig.^(b), where the range is only one 
tenth of the dot size, the splitting between the and 
(— i) classes is still very clear, but the scale of the energy 
levels is much closer to the non-interacting case. Also, 
and more importantly, there is only one singlet-triplet 
crossing in the ground state, as opposed to the three 
crossings of the Coulomb case. For smaller values of oT 1 
these crossings are completely suppressed. 



B. Partition function and susceptibility 

In this subsection we consider the orbital magnetiza- 
tion and magnetic susceptibility of the interacting two- 
electron system. The partition function Z(B,T) = 
2y| e -/3H(S)| ^ _ \/k B T) can be computed from the 
energy levels. The magnetization m(B) and the mag- 
netic susceptibility x(-B) = dm(B) / dB can be calculated 
from 



m(B) 



1 OF 
A dB 



1 dlog Z{B 7 T) 



[3A 



dB 



(16) 



where F is the Helmholtz free energy and A = L 2 . 

Figure || shows the magnetization m(B,T), the sus- 
ceptibility \{B,T) and the interaction-induced suscep- 
tibility x mt — X ~ X^ as a function of B for a = 
(Coulomb potential) and a -1 = L/10 (short-range inter- 
action), is t ne susceptibility of the non- interacting 
system. Both x(_B,T) and x mt (B,T) are expressed in 
units of the Landau susceptibility \xl\ — e 2 /\2irm*c 2 
and the temperature is expressed in units of the mean 
level spacing A. 

On the average, x(B, T) is diamagnetic, as in the non- 
interacting case (see x on the inset). However, the 
exchange-induced singlet-triplet crossings of the ground 
state level contributes paramagnetic fluctuations of or- 
der of ~ 3xl at low temperatures. As discussed in the 
previous subsection, these crossings tend to disappear for 



short range interactions. For aL = 10 only one crossing 
survives. 

For very low temperatures, only the ground state con- 
tributes significantly to the partition function and F « 
Eq. In the region close to the crossing, dEo/dB is dis- 
continuous with a negative curvature, which explains the 
paramagnetic peak in x{B)- For higher temperatures, 
the "anti-crossing" (positive curvature) of the first ex- 
cited state tends to compensate the ground state crossing 
and the fluctuations are attenuated. 

Figure ^ shows the behavior of the interaction-induced 
susceptibility at zero magnetic field x mt (B — 0,T) for 
the Coulomb the and short-range (a -1 = L/10) po- 
tentials and different values of the strength parameter 
Vq. The results show paramagnetic and diamagnetic 
phases as a T increases, with a diamagnetic minimum 
at T w 2A. For low values of Vq the susceptibility is 
again paramagnetic for high T. This type of behavior has 
been observed Jaefore for weak contact (Dirac delta) type 
of interactions!! Our results show that the existence of 
paramagnetic and diamagnetic phases in x mt (B — 0, T) 
might be more general, and not so sensitive to the type 
of interaction. 



V. DISCUSSION 

This work was motivated in part by the experimen- 
tal results of Levy et alQ, who measured the magnetic 
properties of an ensemble of square dots containing a 
few thousands electrons each in the balistic regime. Our 
objective in this paper was to understand the effects of 
the residual electronic interaction in the simplest possible 
case, i.e., that of only two electrons. We simulated the 
shielding of the bare Coulomb force by using an exponen- 
tial type of cut-off, like that of the Yukawa potential. The 
range of the interaction, a, was considered as a free pa- 
rameter. The size of the dot, which controls the relative 
strength of the potential, Vq, with respect to the kinetic 
energy, acts as a second parameter. Our results show that 
both Vo and a are important to determine the proper- 
ties of the energy spectrum and of the probability pro- 
file of the ground state. We showed, in particular, that 
short range interactions may suppress the appearance of 
Wigner molecule type of states even for strong interac- 
tions (Vo large). The range of the interaction also affects 
the magnetic susceptibility of the system. Short range in- 
teractions might inhibit singlet-triplet oscillations of the 
ground state, suppressing the paramagnetic fluctuations 
of x(B)- Finally we have shown that the part of the 
susceptibility induced by the interaction presents para- 
magnetic and diamagnetic phases as a function of the 
temperature, in agreement with the results obtained by 
a pure Dirac delta interaction!! 

It is difficult speculate at this point if our results point 
to an explanation of the slow decay of the susceptibil- 
ity with the temperature, as observed experimentally 
by Levy et al. In a system with many electrons, the 
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ground state energy oscillates as a function of B even 
without electronic interaction. This type of oscillation 
is only due to the boundary, and is responsible for the 
paramagnetic susceptibility of the gas. Our results indi- 
cate that, besides this "boundary-induced" effects, the 
electronic interaction is responsible for further oscilla- 
tions, this time between singlet and triplet states. These 
"interaction-induced" oscillations also contribute to the 
susceptibility. In the present case of two electrons this 
contribution turns out to be larger than that of the non- 
interacting case. Besides, the curves in Fig. ^ show that 
X mt (B — 0, T) is sensitive to both Vq and a. The peak of 
x lnt (B = 0, T) at T fa A leads to a slightly slower decay 
of the overall susceptibility. Although this result might 
be peculiar of few-particle systems, there are similari- 
ties between our findings and those obtained from semi- 
classical analysis .and RPA perturbation theory in the 
high-density limitcl, such as the diamagnetic minima in 
X mt (B = 0,T). This might indicate that the interaction 
indeed plays an important role in behavior of the sus- 
ceptibility with the temperature, although calculations 
with more electrons should be conducted to confirm this 
conjecture. 
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APPENDIX A: CALCULATION OF THE 
INTEGRALS 

In order to calculate the matrix element ( |l4| ) we first 
notice that the integrand can be decoupled using rela- 
tive and center-of-mass coordinates. Since the masses 
are equal (we set mi — mi = 1 for simplicity), we have 
f = fi — f\ and R = (R 2 + Ri)/2, When the sine func- 
tions are written as a sum of complex exponentials 



s'm2nurxi = y~] _j 



sin 27171-?/, = J2 m =-i~ 2l 
the integrand in Eq.(H4) becomes 



2i 
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e i( a i+a 2 )X e i(p 1 +02)Y e i(a 2 -ai)% gi(ft-ft)f ^_ 



where 



„ 7T7711 ,. 



OL 2 



L e 2 



Pi — —Vi + —Vi j Pi — —Vi + —V2 



(Al) 



(A2) 



The integrals on the Center-of-Mass coordinates can 
then be readly evaluated, at the cost of working on the 
complex plane. However, one should note that the limits 
of integration of (x,y) and (X, Y) are not independent. 
In the rest of this development, we will take L = 1 for 
the sake of simplicity. This is equivalent to perform the 
calculations on a adimensional variable -jr. In order to 
consider specific sizes for the dot, one has to include an 
additional factor of L multiplying the whole integral and 
scale a accordingly. 

The change of variables (x\,y±,X2,y2) — > (x, y,X,Y) 
leads us to four different sets of integration limits (the 
quadrants on the (x,y) plane). Next, we show the cal- 
culation of the integral on the first quadrant (7i). The 
calculation on the other quadrants (12,13,14) is analo- 
gous. The total integral is then I = I\ + 12 + 13 + Ia- On 
the first quadrant we have 



i 1 = y2 — " ' e ' 2 ^ m ' ' ' — x 



256 

y e -a^x 2 +y 2 



\/x 2 + y 2 J 
(A3) 

where 

dAi = dx dy J * dX J * dY . (A4) 
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After doing the integrals on the center-of-mass coordi- 
nates we are left with 



(ei---e 2 )(t?i---??2) 
256 



f{x,V)—7 = da;*/ 



(A5) 

where the exact format of f{x,y) depends on whether 
we have (a\ + aa) or (ft -I- ft) equal to zero or not. We 
are going to write it down explicitly in a moment. 
A transformation to relative polar coordinates yields: 



(ei---e 2 )(??i---7?2) 
256 



ri(6) 



f{r,6)e~ ar drdO 

(A6) 



J 



where r\{6) = l/cos6> for < < f and n(0) = 
1 / sin for j < 9 < J. The integral over r can be done 
analytically and we are left with a set of integrals over 9: 



(ei---e 2 )fa---??2) 
256 



(A7) 



As one can see from eq. (A3), the explicit form of F(6) 
depends whether (a\ + 0:2) and/or (ft + ft>) are equal to 
zero or not. We now present the specific form of F(9) for 
all different possibilities. 



1. (ai+a 2 ),(ft +ft)^0 



F(0) 



(-1) 



(ai + a 2 ){0i + ft) 



Z o (0) 



Z d (0) 



(A8) 



2. (ai + a 2 ) = (ft + ft) = 



*W = - (C ° S ^ ing) [(^ (^ -!) + !)] + [(e™ - 2*n + 2) - 2)] (A9) 



Za{0) 



zl{6) " ' " 4(0) 

3. (ai + a 2 ) = and (ft + ft) / 



z b {9) 



z c {0) 



cos 9 



o{0 



i+M (e^(z bri -l) + l) _ (e^(^ri-l) + l) 



(A10) 



4. (qi +q 2 ) / and (ft + ft) = 



F(9) c fa + ^ (^ ri -l) _ («- ri - 1) 



Jq 1+ « 2 ) ( e * dri (^ri-l) + l) _ (e^(z e n-l) + l) 



r 



(AH) 



where 



Z a (8) = 


—a — 


i (a% cos 6 + ft sin 9) 


Z b (9) = 


-a + 


i [a% cos 9 + ft sin 9) 


Z c (9) = 


—a — 


i (ot\ cos 9 — 02 sin9) 


Z d {6) = 


-a + 


i («2 cos 9 — ft sin 9) 


Za{8) = 


-a + 


%{a.2 cos 9 + ft sin 9) 




-a + 


iia.2 cos 9 — 01 sin 9) 




-a + 


%{a.2 cos 9 + 02 sin 9) 


z d (9) = 


-a + 


i(02 sva.9 — a± cos 9) 


Ze(9) = 


-a + 


i (02 sin 9 + 0*2 cos 9) 



All non-trivial integrals over 9 have been performed 
numerically. 
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FIG. 1: Energy levels as a function of Vb for different values of the reach parameter a. (a) a — (Coulomb interaction), (b) 
aL — 1 and (d) aL = 10. Inset: avoided crossings on the (+l)-singlet (solid line) and (— l)-singlet (filled squares) classes 




FIG. 2: Electronic density for different values of a 1 and L. For each line, the reach parameter a 1 is fixed for increasing 
values of L (from left to right L — 10, 100 and lOOOnm). From top to bottom, we have a" 1 = 10, 100 and lOOOnm. The Wigner 
Molecule state is recovered for a -1 = L = \\xm 
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FIG. 3: Top: Electronic density (right) and probability density with one of the coordinates on a Wigner molecule peak 
(|*^ C) (ri = (0.2, 0.2), r 2 )| 2 ) (middle) showing a singlet-like spatial correlation (right). Bottom: Same for an excited triplet 
state [(-fi)-Class] 




(b) 

{^'"I'iii"' ■■wii C • 

:n;:.-'**«S: y**„„w 

1 1.5 k; ; -s-' - -» ":•;.;::.•■;.< 

ir ; x 1 '"' •0* ,l "s«: !,,! "'.ss»' 



3:Ti 



(C) 

;.;i,.,i:;;,.,„«K::":""::ii 

>><„ >C-.. 

■.^jiiM"""** ■ Class ■ ! isrrjlti; 

• Class +1 (triplet) 
'Class +i (singlet) 
' Class +i (triplet) 
) Class -1 (singlet) 

• Class -1 (triplet) 
■-■Class -i (singlet) 



BIT) 



FIG. 4: Energy levels as a function of B for L = 200nm and different values of the reach parameter a 1 . (a) non-interacting 
(Vo = 0), (b) a" 1 = L/10 and (c) Coulomb potential (a = 0) 
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FIG. 5: Magnetization m(B), magnetic susceptibility x{B) an d interaction induced susceptibility \ m = X~X as a function 
on the magnetic field B for a dot size of L = 200nm and different values of the potential reach parameter a (Top: a = 
[Coulomb], Bottom: a -1 = 20nm). For low temperatures (T < A/2), the system features interaction-induced paramagnetic 
peaks. Inset: The non-interacting magnetization and susceptibility. 
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FIG. 6: Interaction-induced zero field susceptibility Xjj=o(T) as a function of the temperature for a Coulomb potential (left) 
and a Yukawa potential with a -1 = L/10 (dir.) and different values of the potential strength Vo- 



